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Periodic Mean-Field Solutions and the Spectra of Discrete Bosonic Fields: 
Trace Formula for Bose-Hubbard Models 
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Institut fiir Theoretische Physik, Universitdt Regensburg, D-93040 Regensburg, Germany 

We consider the many-body spectra of interacting bosonic quantum fields on a lattice in the 
semiclassical limit of large particle number N. We show that the many-body density of states can 
be expressed as a coherent sum over oscillating long-wavelength contributions given by periodic, non- 
perturbative solutions of the, typically non-linear, wave equation of the classical (mean-field) limit. 

To this end we construct the semiclassical approximation for both the smooth and oscillatory part 
of the many-body density of states in terms of a trace formula starting from the exact path integral 
form of the propagator between many-body quadrature states. We therefore avoid the use of a 
complexified classical limit characteristic of the coherent state representation. While quantum effects 
like vacuum fluctuations and gauge invariance are exactly accounted for, our semiclassical approach 
captures quantum interference and therefore is valid well beyond the Ehrenfest time where naive 
quantum-classical correspondence breaks down. Remarkably, due to a special feature of harmonic 
systems with incommesurable frequencies, our formulas are generically valid also in the free-held 
case of non-interacting bosons. 


I. INTRODUCTION 

The full quantum mechanical solution of the problem 
of interacting particles gets exceedingly complicated with 
increasing particle number, and even for a generic single¬ 
particle problem in the limit of large excitations. Hence 
there has been the quest for devising versions of the quan¬ 
tum formalism where classical input can be used to pre¬ 
dict the outcomes of observations keeping intact concepts 
like superposition of states and summing amplitudes in¬ 
stead of probabilities as embodied in the kinematical 
structure of quantum mechanics. A natural benchmark 
for the use of classical objects in quantum mechanics is 
the ubiquitous and defining presence of interference phe¬ 
nomena in the quantum world. 

One attempt to search for quantum effects using only 
classical information consists of following the time evolu¬ 
tion of quasiclassical, coherent quantum states with the 
sharpest distribution of momentum and position allowed 
by quantum mechanics. The time evolution of minimal 
wavepackets is then approximated for short times by a 
rigid motion along the unique classical trajectory fixed 
by the initial expectation values of position and momen¬ 
tum ((7(0)), (p(0)). For times shorter than a usually short 
characteristic quantum scale, the Ehrenfest time, expec¬ 
tation values are given simply by the classical values 
q{t),p{t) determined by the unique solution of the clas¬ 
sical equations of motion. This approach breaks down 
when the time evolved expectation values of {q{t)), {p{t)) 
are insufficient to recover even approximately the time 
evolved wavepacket. This happens when different sectors 
of the originally well localized wavepacket start superim¬ 
posing with each other and produce interference patterns. 

The failure of this approach lies in its direct use of 
classical concepts, as there is no simple way to modify 
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classical mechanics in order to account for interference 
phenomena, particularly if the approach is fixing a unique 
classical trajectory. The early realization that quantum 
phenomena can indeed be explained in terms of inter¬ 
fering amplitudes between classical paths instead of the 
classical trajectories themselves marks the beginning of 
the semiclassical program (see [l| for a historical review). 
It took, however, almost fifty years until Gutzwiller pro¬ 
vided a complete and rigorous derivation of the semiclas¬ 
sical approximation to the quantum mechanical propaga¬ 
tor, the starting point of modern semiclassical methods 

3. 

For first-quantized, single-particle systems, 
Gutzwiller’s result for the quantum mechanical density 
of states (DOS) 

p{E)=Y,S{E-E^), ( 1 ) 

n 

where E^ are the eigenvalues of the Hamiltonian, has the 
generic form 

p{E) ^ -p{E) + p{E) (2) 

in the formal limit h ^ 0. Here, the smooth part p{E) is 
purely classical, in that it is related with the phase space 
volume of the classical energy shell, also known as the 
Weyl term. Remarkably, quantum fluctuations responsi¬ 
ble for the oscillatory part p{E) are also given in terms 
of classical quantities, though encoded in a subtle way in 
the periodic solutions of the classical equations of motion 
of the corresponding classical system. A key distinction 
between these two contributions to the semiclassical den¬ 
sity of states is that while p{E) is analytical in h and 
therefore admits a power expansion around h = 0, this 
is not the case for the oscillatory contribution, namely, 
p{E) cannot be approximated by any finite-order expan¬ 
sion in h. 

It is natural to ask which modifications, technical or 
conceptual, are required to take the semiclassical pro- 
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gram into the realm of many-body systems where corre¬ 
lation effects due to both interactions and indistinguisha- 
bility additionally appear. Here we face a unique as¬ 
pect of many-body systems, namely, that one can choose 
between two, equivalent but conceptually quite different 
approaches. On the one hand, one can generalize first- 
quantized techniques to many-body systems by extend¬ 
ing the number of degrees of freedom and using projector 
techniques to select the states with the appropriate sym¬ 
metry under permutation. The associated semiclassical 
approach is then based on interfering classical paths in a 
multidimensional space, supplemented by boundary con¬ 
ditions and/or extra classical paths joining initial classi¬ 
cal configurations with the ones obtained under permuta¬ 
tions of the particle labels. This first-quantized approach 
has been successfully applied to the helium atom, a pro¬ 
totype of a strongly interacting few-particle problem Q , 
and has been used so far to derive formal results for the 
symmetry-projected trace formula [ 1 , @|. Moreover, the 
simplifying assumption of a unique mean-field potential 
fixing the classical dynamics and making it essentially 
non-interacting has lead to the notion of shell effects 
and its semiclassical interpretation 0, @1 . Furthermore, 
a complementary semiclassical approach to the smooth 
part of the DOS in many-body systems has shown to be 
surprisingly accurate Q and allows for going beyond the 
independent particle model. 

The above mentioned approaches have in common that 
the semiclassical limit of high energies (/i —?> 0) is taken 
for fixed, though possibly large total particle number N. 
The second option to describe quantum many-body sys¬ 
tems of interacting identical particles is the use of quan¬ 
tum fields. In this approach, the quantum dynamics for 
iV —>■ oo has as classical limit a, typically non-linear, wave 
equation. This new feature is, however, compensated by 
the fact that with the notion of quantum fields, indistin- 
guishability is included by construction into the kinemat¬ 
ics of the state space, instead of by applying projectors 
as in the first-quantized approach. The use of quantum 
fields to describe systems of interacting, identical parti¬ 
cles has another key consequence. It has been shown (and 
it will be apparent from the calculations presented here) 
that for systems with fixed, finite-dimensional Hilbert 
space, the classical limit is actually equivalent to N ^ oo. 
It is then reassuring that, as shown in [l3,Eli, the classi¬ 
cal limit of the second-quantized theory turns out to be 
the mean-field description, which is expected to correctly 
represent the dynamics in the thermodynamic limit. 

In this paper, we follow Gutzwiller’s program and de¬ 
rive rigorously a formula providing the DOS of a second- 
quantized system for large N where the classical limit in¬ 
volves discrete field equations, the mean-field equations 
of the associated discrete quantum field. We construct 
both the smooth and oscillatory contributions to the 
quantum many-body DOS starting with the semiclassical 
approximation to the exact Feynman propagator for the 
quantum field. We show that the many-body DOS arises 
from interference of, in principle, infinitely many, peri¬ 


odic solutions of the corresponding classical mean field 
equations, in close analogy to the periodic orbit contri¬ 
bution to the single-particle DOS. In the case where in¬ 
teractions of strength g are present, our derivation relies 
on the existence of chaotic behavior in the classical limit 
(the ubiquitous presence of chaotic regions in the phase 
space of discrete mean field equations has been addressed 
in 0 )- For this, a suitable scaling of the interaction 
strength g ^ N~^ is usually assumed that keeps its con¬ 
tribution to the total energy U ^ gN'^ proportional to N 
and therefore comparable with the single-particle (hop¬ 
ping) term. In the non-interacting case our derivation is 
valid for systems where the single-particle spectrum used 
to construct the many-body state space is such that the 
energies are non-commensurable, as it is generically the 
case. Our work for the free, non-interacting case opens 
a road to study in a systematic way the combined limit 
N —>■ oo,h^ 0 in infinite-dimensional systems where 
new kinds of classical structures might become relevant 


A central aspect of Gutzwiller’s method [1] is a clever 
choice of the representation where the semiclassical prop¬ 
agator, the key object representing quantum evolution in 
terms of solutions of the classical equations of motion, 
appears as a sum of oscillatory terms given by real ac¬ 
tions. While in the first-quantized scenario this choice is 
naturally given by the position representation, this im¬ 
portant aspect of the semiclassical program has not been 
addressed in the context of quantum fields, where the 
usual choice for constructing the path integral is the co¬ 
herent state representation for which the actions entering 
the semiclassical propagator are complex 0- Thus, here 
we use a different approach. 


The paper is organized as follows. A key point in our 
approach is to generalize the concept of position eigen¬ 
states into the realm of quantum fields, and for this rea¬ 
son in Sec. mi we briefly introduce these objects, together 
with basic definitions of Fock space and creation and an¬ 
nihilation operators. The main technical part of the pre¬ 
sentation, Sec. mil is dedicated to the rigorous deriva¬ 
tion of the many-body DOS for Bose-Hubbard systems 
starting from its quantum mechanical definition and the 
semiclassical approximation of the quantum mechanical 
propagator. Specifically, the smooth (Weyl) contribution 
to the DOS is derived in Sec. IHI Al with its final form 
given in Eq. m, while the oscillatory, Gutzwiller-like 
contribution is derived in Sec. IHI HI and summarized in 
Eqns. (TOD . The all important derivation of the correct 
phases accounting for the existence of focal points, the 
Maslov indeces, is outlined in Sec. IHI Gl We finish this 
technical part with Sec. IHI PI where we explain why the 
extended phase-space approach used here to account for 
the peculiar consequences of number conserving Hamil¬ 
tonians has advantages over a reduced-phase space ap¬ 
proach, although the later is more frequently used to deal 
with continuous symmetries in first-quantized systems. 
As a first application of our formalism, in Sec. II VI we de¬ 
rive a semiclassical expression in terms of intefering peri- 
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odic mean-field solutions that gives the exact form of the 
many-body DOS for a non-interacting bosonic field with 
incommensurable single-particle energies. The connec¬ 
tion of our approach with some existing semiclassically 
inspired ideas and methods in the context of interacting 
bosonic systems is presented in Sec. lYl together with pos¬ 
sible applications. Finally, in Sec. IVIl we conclude with a 
summary of our findings and some of its consequences. 


ni,... ,nL- These states satisfy a^a; |n) = n; |n) and 
1 


|n)=- 


Iflnil 

1^1 


(4) ‘-(at) ‘lo), 


(4) 


More important for the derivation of the trace formula, 
however, are so-called quadrature ei gens tates |q) and |p) 
defined by the eigenvalue equations [1^ 


II. QUANTUM MECHANICAL BACKGROUND 


We will restrict ourselves to quantum fields described 
by a general Bose-Hubbard Hamiltonian with two-body 
interactions. 


H = 


L 

E 

1 ,^2 = 1 


H, 




1 


L 


T 2 / . 

L ^ 2,^3 —1 




(3) 


Here, H = ^ is the hermitian matrix de¬ 

scribing the single-particle motion and the four-fold sum 
describes two-body interactions. Moreover, ai and a| 
are the annihilation and creation operators for the Tth 
single-particle state (or site) satisfying the usual bosonic 


commutation relations 


ai,al 


5w. 


At intermediate steps, we will make use of Fock states 
|n) determined by the (integer) occupation numbers 

I 


- (a; +djt) |q) = qi |q) , 

( ; ; (5) 

which satisfy 

and 

{v\^) =W_—/====HnA'/2piY ( 7 ) 

where denotes the n-th Hermite polynomial. We have 
also the overlap 


L 

(q|p)=n 

1=1 

and the closure relations 


Q^ipiqi 



i= d^q |q) (q| = / d^p |p) (p| . 


( 8 ) 

( 9 ) 


For the derivation of the smooth part of the DOS, we 
will later make use of the asynmtotic formula of the Her¬ 
mite polynomials for large n [8|, 


(q|n) 



( 10 ) 


III. DERIVATION OF THE MANY-BODY 
DENSITY OF STATES 


with 3 denoting the imaginary part and with the trace 
in terms of a sum over Fock states 


The DOS Pn{E) for fixed number of particles N is 
given by the imaginary part of the trace of the Green 
function G{E) over the subspace of the full Hilbert space 
obtained by fixing N = It is given by 

Piv (£^) =--lim %Ar(£’-l-iT?), (11) 

TT rj-^0 


9n{E) 


TtnG (E) = 


} Sr. (n 

“ 1 = 1 


G{E) 



( 12 ) 

Semiclassically, the single-particle DOS is typically split 
up into a smooth part, which stems from short trajecto¬ 
ries and an oscillatory part determined by periodic orbits, 
see Eq. ©• Correspondingly, for the many-body case we 
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will now first derive the smooth part pj^{E). 


with the Weyl symbols of an operator O being defined 
by M 


A. Smooth part 


To this end we first rewrite the sum over all possible 
occupations and the Kronecker delta in Eq. dm) by a 
sum over those occupations, which have the correct total 
number of particles and insert the definition of the Green 
function as a Laplace transform of the propagator K{t) 
in Fock space, 


6 


Weyl 


(q, P) = 




6 


q 



(15) 



gN{E)= 


n 


1=1 


= -/dte^®* ^ (n 


(13) 


n 


n: J2 ni=N 

1=1 


The smooth part of the DOS stems from short paths, 
i.e. from the short time contribution to the integral. In 
order to compute this contribution, we will first evaluate 
the trace and then perform the integration. To this end 
we rewrite the diagonal matrix elements as 


n 


K (t) n\ = TiK{t) |n) (n| = 


(14) 


d^q J d^p ^A:(t)J^^^^(q,p) [|n) (n|]^^yi (q, p) 


Next, we use the asymptotic formula (I10|) for the Her- 
mite polynomials for large n and rewrite the cosines as 
exponentials, yielding four terms, where for two of them 
the exponents from the cosines have the same sign, while 
for the remaining two these signs are different. However, 
in Eq. m replacing Q by — Q is the same as complex 
conjugation, which is only true, if the two signs of the ex¬ 
ponential are opposite. Therefore, the terms with both 
signs being the same have to cancel when performing the 
integral. Thus, the resulting exponent is antisymmetric 
in Qi. Expanding it up to second order in Qi and ne¬ 
glecting the dependence of the prefactor on Qi yields 


^ / dQzexp -2iQi (pi + si^{ni + 5 ) - gf 

]n) (n|]weyl (q-p) ~ n - 

1 = 1 Sl=±l 


5 (pi + Si {ni + ^) - qf'j 
sj{ni + i) - qf 


-HE, 

1=1 Si=±l 


1=1 


I- 

For the Weyl symbol of the propagator, one can use Here, 
the usual short time asymptotic form 


(p,q) = 


H 


(p|q) 


(16) 


(18) 


K(t) (q,p)Riexp 

J Weyl 




E f “ b E 

ll.l2 = l V ^ ^3 = 1 / ^ ^ 

hhhU ^ ^ ^ ^ 

^ ^ is the mean field Hamiltonian corresponding to 
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the full quantum Hamiltonian (jS)). It can be obtained by the simple replacement rule [isini 

a]au ^ ipiipv - ^Sw, ( 19 ) 

with = qi + ipi. 

Inserting Eqns. (HH) and dm) into Eq. da as well as 
replacing the sum over occupations in Eq. da by an 
integral then yields for the smooth part of the resolvent 


gN{E) = 




dt d^q d^pexp 


h . 


E- (p,q) 


t U q 2 + p 2 - iv - 


0 


and thus for the smooth part of the many-body DOS 

L 


-PN{E) = J d^q J dV {e - (p, q)) S - N - 


( 20 ) 


( 21 ) 


As in the single-particle case, the smooth part is given by 
the phase space volume of the fV-particle energy shell. 


B. Oscillatory part 


To compute the oscillatory part g^iE), we start again 
from Eq. (IT^ and rewrite the Kronecker delta, to get the 
resolvent as 


9n{E) 


277 / 

/daE(n 

— ia E ^ 

e V 7^1 ' )g{E) 

0 ^ 




( 22 ) 

The oscillatory part of the DOS, which we are in¬ 
terested in here, can be obtained from a semiclassical 
approximation of the Green function by computing the 
trace using a stationary phase approximation. However, 
for Eock space, the stationary phase approximation is not 
applicable, since the trace in Eock states is given by a sum 
rather than an integral. On the other hand, in [l^, [ll| a 
possible way to circumvent this problem has been shown. 
This is by again using the quadrature eigenstates |q), |p). 
Inserting them to the left and to the right of the Green 
function yields 


gN {E) = 

2tt cxd oo 


( 23 ) 




da j d^q j d^p'^ip G (E) qyqln) 
-ia j 


0 OO OO 


X ( n 


p • 


Using the completeness relation of the Hermite polyno¬ 


mials then yields 


E(qi“) 



A exp [2qipi + {qf +p‘f) sing] - if} 
-v/tt cos a 

7—1 V 


( 24 ) 


Next, one has to find an expression for the Green func¬ 
tion, which is related to the propagator by means of a 
Laplace transform. 


G'(p,q;E) = G (E) qj = — J dte'>^^*K {p,q,t). 

0 

r-i . . . . 

In [ll| , a semiclassical approximation for the propagator 

has been found, which is given by 


K (p,q,t) = E 


det ^ 


2ttH dpdq 


en 




( 26 ) 


where the sum runs over all mean-field trajectories (non¬ 
linear waves) 7 given by the solutions of the equations of 
motion 




(■ 0 *(i)^ ^^(i)) 


( 27 ) 


and the boundary conditions 


5 Ri/j( 0 ) = q, 

3"0(7) = P- 


( 28 ) 


Moreover, the phase each trajectory contributes with is 
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given by its action 


bians, with the period of 7 , 




R-y= dt' 


■ 0 (t')) 


— 2p • 


(29) 


and the Morse index 

For later reference, we state the derivatives of the ac¬ 
tion with respect to p, q and t\ 

2n^i){t), (30) 

H(mf) (^*(o),,/,( 0 )) = 

In order to determine the oscillatory part of the many- 
body DOS, the time integration in Eq. (I25|) can be eval¬ 
uated using a stationary phase approximation. The sta- 
tionarity condition then selects those trajectories which 
have energy E, 


dRy 

dp 

dRy 

dq 

dRy 


dt 


[Ry -|- Et] = E — Ey = 0 


(31) 


d{^i,{0),Ty) 

dip,E) 


— det 


= det 


/ dj^^m^d{p^\ 

V dip,Ty) d(p,E)J 
/dQtpiO)\ dTy 

V c>p J dE' 


With this, the semiclassical Green function Eq. ((25)l is 
given by 


G(p,q;E) = 
1 1 


(33) 




\/2'Kh 


-T. 


\ 


a'^w^ 


det 

I a^w-y a^Wy 
\ aEdp aE'^ 




jrWy-itly] 


with ^y = fly + sign{dE/dty)/2 and 


Wy=Ry+ ETy =2h 5i/j(t) • 5it/;(t)dt - 2ftp • ^ip{Ty) 


satisfying 


(34) 


dWy 

dp 

dWy 

dq 

dWy 

dE 


= - 2hUip{0), 
= - 2hSsxp{Ty), 
= Ty. 


(35) 


In order to compute the semiclassical prefactor of the 
Green function, one can use the standard trick for Jaco- 


Thus, in the semiclassical limit the oscillatory contribu¬ 
tion of the resolvent ( 1 ^ is given by 


qn^E) = 


27ri?i 


dae 


— \OlN 


<17 / <i 7 E 




det 


a'^Wy a'^Wy 

aqSp 

a'^Wy 


aEdv 


dcidE 

a'^Wy 

aE'-‘ 


h ^ 2 ^ 2 




i[^2q-p+(^q^4-p^) £ 


V2TrH yjE cos ot 


(36) 


The integrations over p and q as well as a will again be The corresponding stationary phase conditions for the 
performed in stationary phase approximation. integrations over p and q read 


J 




[2q • p -I- (q^ -I- p^) sin a I 
cos a 

[2q • p -I- (q^ -I- p^) sin a I 
cos a 


= -2 




= -2 


^^{Ty) - 


qsma 
cos a 


p sm a 
cos a 


cos a 


cos a 


= 0 , 
= 0 , 


(37) 

(38) 


which can be combined into the more compact condition 
t/.(r^) = i/i(0)e-'“. (39) 


Equation (I39|) implies that the resulting trace formula 
will be given by a sum over pseudo-periodic orbits, [l^ 
for which the associated classical nonlinear waves after a 
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certain pseudo-period differ from their initial values 
by a global phase a. 

Representing the classical nonlinear wave solution in 
terms of its amplitude and phase, 


= \/ni 


(40) 


the resulting stationary phase solution for the a- 
integration in Eq. (l36l) is given by 


S^ = W^ + 


[2q • p + (q2 + p2) si 


Sin a 


cos a 


= hj e{t) ■ h{t)dt + hn{Q) ■ [0(0) - 0(T^)]. 


(41) 


In the last term, which originates from a partial integra¬ 
tion, n(T-y) = n(0) has been used. Its dependence on a 
is determined by 


n(0) • [0(0) - 0(T^)] = N^a + 2^ ^ ni{0)ki, (42) 


1=1 


with 


2L — 1 of the stationary phase conditions, Eqns. (I37I38I) 
suffice to satisfy all of them. Due to the Noether theo¬ 
rem, there is a continuous symmetry for each conserved 
quantity. Here, this continuous symmetry is given by 
the ?7(1) gauge symmetry, i.e. the freedom to chose an 
arbitrary time-independent global phase 9. 

Moreover, as in Gutzwiller’s original derivation Q of 
the single-particle trace formula, the starting point of the 
pseudo-periodic orbit can be chosen at any point along 
the orbit. Thus, there remain two integrations, which are 
the integrations over all pseudo-periodic orbits belonging 
to the same continuous family of trajectories, that have 
to be performed exactly. 

Eor single-particle systems, the trace formula for 
chaotic systems with additional continuous symmetries 
has been studied in [l^. The evaluation of the semi- 
classical prefactor of the trace formula presented there 
can, to a large extend, be carried over straight forwardly 
with minor modifications in order to correctly account 
for the fact that the orbits in the case studied here are 
not strictly periodic. Therefore, here we will only show 
the steps, which have to be altered and refer the reader 
to Ref. [13 for more details. 






(43) 


the (time independent) number of particles defined by 
the trajectory and fci,..., /cl being integers. 

Taking a closer look, one recognizes that due to the 
conservation of the total number of particles already 


After transforming in Eq. (1361) the integration vari¬ 
ables q,p locally to qy, q_L, P||, P_l, where the (two- 
dimensional) parallel components run along directions of 
the continuous families, i.e. along the trajectory and the 
direction of the global phase 9, while the remaining ones 
are perpendicular to these, the integrations over the per¬ 
pendicular components as well as pii yield 


2 ir 


i«(E) = - 


da / d^ciii Yi. : 


:,-i“(^+T) 




cos a 


3(ad.(o),T,) 

a(p,B) 


r 


(0) sin a. d'isipj_{0) 1 

Oq_L cos a dp cos a 

_ 1 _ sin a. 

Oq_L cos a dp cos a 


-1 


:» h ^'y ^ 2 ) 


(44) 


where = {N+ — NJ)/2 is the difference between 
the number of positive and negative eigenvalues of the 
(2L — 2) X (2L — 2) dimensional matrix appearing in 
the semiclassical prefactor. Note that the sum runs over 
pseudo-periodic orbits with the initial global phase and 
the initial position within the orbit chosen by the integra¬ 
tion values. Alternatively, one can also refer to the sum 

I 


over 7 as a sum over families of pseudo-periodic orbits, 
where one is free to choose the initial global phase of the 
reference orbit, which is used to compute its contribution. 

Leaving the calculation of the determinant in Eq. (HU) 
to Appendix|3 the trace of the semiclassical Greens func¬ 
tion is given by 


2 ir 


~9n{E) = - 


da / d^^ii ^ ■ 


-ia(jV-Sl) 




det (M.y — I) 


ae 

dN.y 


dift 


(^11 (7^7) e‘“) 


5(r^,0) 




(/X7-+I^7+L77(a)) 


(45) 


Here 9 is the initial phase of the trajectory, 

f I if 5 < a < 

0 else 


Stt 
2 ’ 


and 


(46) 


9(q^ 


5(5R(i/>_L(T^)e-),Pi) 




u( 0 )e-'“)) 
a (5? (i/>^(T,)e‘“),5(^_L(7’7)e'“)) 
(0)),3(tA^ ( 0 ))) 


( 47 ) 




























is the stability matrix for the pseudo-periodic orbit. In 
view of Eqs. (I28I39I) . for the pseudo-periodic orbit 

5 (t/;( 0 )e-‘“) =p, 

/ . (48) 

K = q 


holds. 

Hence, the matrix in Eq. (14711 is indeed the many-body, 
field-theoretic analogue to the monodromy matrix ap¬ 
pearing in the usual Gutzwiller trace formula Q. 

Now, in Eq. (|4^ . the last determinant can be used 
in order to transform the integration over qy into inte¬ 
grations over the propagation time and the global phase. 
Again in view of Gutzwiller’s derivation [3| one has to 
correctly account for repetitions of each primitive pseudo- 
periodic orbit when evaluating these integrals. These 
primitive pseudo-periodic orbits are obtained by find¬ 
ing the largest possible, but finite integer m > 1 for 
which 'ip{T.y/m) = 'i/j( 0) exp(—ia/m). Then 'ip{t) obvi¬ 
ously still satisfies Eq. (IMl) . However, after the pseudo¬ 
period Tppo = T^/m the primitive orbit is repeated but 
with a different global phase. Thus when naively inte¬ 
grating the global phase from 0 to 27r and the time from 
0 to Tly, one and the same orbit is counted m times. 

On the other hand, as discussed in Appendix m for 
a given pseudo-periodic orbit, the primitive phase dif¬ 
ference a, i.e. the phase difference after the primitive 
pseudo-period, is unique. That is that any time T*, 
for which ipiT*) = ■0(0) exp(—ia*), has to satisfy T* = 
iTiTppo with m G N. Thus obviously also a* = ma. 

Therefore, 


where Tpo is the flying time of the orbit, which may be 
any integer multiple of the primitive period Tppo, defined 
as the smallest time, for which Eq. (15^ is satisfied, and 
Qfpo is an arbitrary global phase depending only on the 
trajectory. 

The argument of the cosine is given by the classical 
action 


-LpO 

Spo{E) = h J 9 {t) ■ n (t) dt + 2Tihkpa 


(53) 


and the (integer) Maslov index 


T / \ 1 • ^-^7 

^po — Mpo “1“ ^po (oipo) ^ 


(54) 


We would like to give a final remark about the appear¬ 
ance of the global phase difference Upo- When consider¬ 
ing the orbits in the reduced space, where not only the 
number of particles is fixed but also the global phase is 
set constant, they would be strictly periodic. However, 
as it was already remarked in [T^ . a trajectory, which is 
periodic in reduced space may not be periodic in the full 
space. 

On the other hand, one might have expected this be¬ 
havior already in advance, since even if the nonlinear 
wave at final time differs from the initial one by a global 
phase factor, the following time evolution is again the 
same as the initial one. 


C. The Maslov index 


f r|2 

(95R (^011 (^)e‘“) 

/ ^ 'll! 

d{T„e) 


2 'kTj 

m 


SttTppo. (49) 


The last remaining integration over a can straightfor¬ 
wardly be computed in stationary phase approximation. 
The stationarity condition selects those trajectories, for 
which the given number of particles is related to the 
total number of particles, 


While Eq. (IMl) in principle yields the correct Maslov in¬ 
dex, it is not very helpful when calculating it in practice. 
A more useful formula can be obtained by not performing 
both integrations over p and q in Eq. (1551) together but 
one after the other. For instance, if p is integrated out 
first, the intermediate result for the resolvent is given by 

271- 

ffN(E) =-i^ /dae-“(^+*) / d^q 


N^ = N+-, 


(50) 


27 ri?i (— 27 rifi) ^ 


z 


/ 


a^w^ \ 

however, when evaluating the integral, one should keep 


det I 

dqdq' 

a^w^ 

dcidE \ 

a^w^ 1 

in mind that a = 9 — 9{Tj), where 6{Tj) is the global 

7 \ 

\ 

\ dEdci' 

aE^ / 


(55) 


ei^. 


phase at final time. 

Finally, the oscillatory part of the many-body DOS for 
fixed total number of particles then reads 


where now the trajectories satisfy the boundary condi¬ 
tions 


Pn{E) = ^ 


po 


- .cos(^Spo{E) - CTpo^ 


7rhy/\Mpo - 1 | \fi 


(51) 

Here the sum runs over the families of pseudo-periodic 
orbits satisfying 


0(rpo)=0(O)e-‘“-, 


(52) 


3?0(O) = q 
5R0(T^)e'“ = q, 

and their actions are given by 

Wj = W^ + -^ 130(Ty) • q -f [q^ -f (S0 (Tt,)) 

COS Gi ^ L 


(56) 


sina| 

(57) 
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Performing the remaining integrals in stationary phase 
approximation (except of those along the trajectory and 
along the global phase) must finally again yield Eq. (ICTI) . 
However, this way following [2^, the Maslov index is 
given by a sum of two terms 


o-po = Mpo + (58) 

where Upo is increased and decreased by one every time 
the determinant of 


/ 93'i/)j_(t)e'“ 53'i/)j_(t)e‘“\ 


changes sign as a function of t. In fact, is not an 
invariant property of the pseudo-periodic orbit, but de¬ 
pends on the choice of the initial point, also depends 
on this choice and is determined by the zeros of the de¬ 
terminant of ■ This index can be deter¬ 

mined as follows [^: When shifting the initial point 
along the orbit, a caustic, which is the point at which 

( ^*( 0 ) ) zero, can appear or disappear. At such 

a point, is incremented or decremented by one. This 
way, (Tpo is independent of the choice of the initial point. 


D. Reduced vs. Extended phase space approaches 

We conclude the presentation of the derivation of the 
trace formula for second-quantized many-body systems 
with a remark concerning the implementation of the 
gauge symmetry responsible of the conservation of N. At 
first glance one may think that our choice of using the 
periodic orbits in the extended phase-space, thus render¬ 
ing them pseudo-periodic, leads to substantial technical 
complications compared with a construction based on pe¬ 
riodic orbits in the reduced phase-space fixed by the total 
number of particles. The classical mean-field equations 
get, however, extremely involved when one explicitly uses 
the conservation of N to reduce the dimensionality of the 
problem, as it can be easily seen from Eqns. 63]) and 
(HU. Explicit use of N as an external parameter leads 
then to equations of motion which are non-polynomial 
in the fields, thus rendering both analytical and numer¬ 
ical calculations much more difficult already in the non¬ 
interacting case while, with our choice, the simplicity of 
the mean-field equations is preserved. This fact will be 
more evident and crucial in the next section where we 
completely and exactly solve the free-held case, some¬ 
thing which is possible because of the strict linearity of 
the problem when formulated in the language of pseudo- 
periodic orbits. 


IV. THE FREE FIELD 

The trace formula (ISTl) for Bose-Hubbard systems hnds 
its most natural application in the case where the mean- 
held equations display (discrete) held chaos going along 


with isolated unstable periodic solutions. Technically, 
this stems from the essential step where the integrations 
involved in the calculation of the trace are performed in 
stationary phase approximation. Implicitly, we are as¬ 
suming that periodic orbits are isolated and do not come 
in continuous families. The presence of continuous fam¬ 
ilies oEperiodic orbits is a hallmark of classical integra- 
bility [ij , and therefore the trace formula cannot usually 
be applied to integrable systems. 

In order to study the possible application of the trace 
formula m in the non-interacting limit of a discrete 
bosonic held we must check whether this limit corre¬ 
sponds to a classical integrable system or not. In view 
of Eqs. 61 ^n.d (fTOl) . the mean held Hamiltonian corre¬ 
sponding to the quantum mechanical free held Hamilto¬ 
nian, 


H = '£H,,ala„ (60) 


is given by 

^(MF.Free) ” y) (61) 

ij ^ Z 

where the term 1/2 came from the Weyl ordering of oper¬ 
ators implicit in our derivation of the semiclassical prop¬ 
agator. 

First, we will show that ij(MF,Free) gg^ gf 

independent constants of motion, implying by dehnition 
integrability. In exact analogy with the quantum case, 
we consider a transformation 

V'i = XI ( 62 ) 

X 

which is canonical if and only if the matrix u with entries 
Ui^ is unitary. It is a simple exercise to show that if the 
matrix u diagonalizes the matrix H, i.e, 

= e^S^y-i (63) 

b 

then the functions 

% (^*: i’) ■= X ^0 > for X = 1,..., L (64) 

b' 


constitute L independent constants of motion under the 
Hamiltonian how induced by ijl^F.Free)^ These classical 
phase-space functions are the obvious classical analogues 
of the quantum mechanical number operators counting 
excitations in the eigenstates of the single particle Hamil¬ 
tonian. 

Since the free mean-held Hamiltonian is integrable, the 
trace formula should in principle be modihed to account 
for the continuous families of periodic orbits typical of 
integrable systems. Remarkably, it turns out that the 
non-interacting limit of a quantum held theory is not 
typical at all. The reason is that, as it is obvious from 
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the quadratic dependence of ij(MF,Free) on the 

canonical variables ip* and i/j, the free field is not only 
integrable but it is actually harmonic. Harmonic systems 
are not generic integrable systems. In fact, depending on 
the number-theoretical relation between the energies 
of the single particle orbitals, they share some funda¬ 
mental properties of the chaotic case. In particular, if 
the single-particle energies are not commensurable (the 
generic situation for a randomly chosen matrix H), the 
periodic orbits of the system are actually isolated. To 
understand this we focus on the solutions of the classical 
limit, which is just the single-particle, linear Schrddinger 
equation 

= X! (65) 

3 

with solution 

ip{t) = e~^^*ip{Q). (66) 

Note that the eigenvector of H 

(67) 


single-particle problem. For fixed energy they are obvi¬ 
ously discrete, and therefore isolated. 

It is important to stress that, being simply the classical 
limit of the theory, there is no physical reason whatso¬ 
ever to prefer a normalized solution of the equations of 
motion. In fact, each eigenvector defines a com¬ 
plete, continuous family of periodic orbits with norms 
that vary continuously. As it will be clear below, this 
continuous family corresponds to the expected continu¬ 
ous dependence of the action with the energy E. For 
fixed E in the trace formula, a specific value of the norm 
I will be selected. 

Note here that in this special case of zero interactions 
the integration along the periodic orbit is actually the 
same as the one along the initial global phase, due to 
Eq. (IMI), namely the time evolution for a periodic orbit 
is simply a change in the global phase. Thus, contrary to 
Sec. lIII B1 only one integration, namely the one along the 
trajectory, has to be performed exactly. The evaluation 
of the integrals in Eq. (l36l) is then strictly equivalent 
to the standard derivation of the trace formula Q and 
yields for the oscillatory part of the many-body DOS for 
non-interacting systems 


with eigenvalue 


( 68 ) 

defines a family of periodic orbits with fundamental fre¬ 
quency 

v(x)(t) (69) 

To show that these are the only periodic orbits of the 
system and that they are indeed isolated we note that, 
because of linearity, Eq. (I5S1) can be expressed as a linear 
combination. 


t/;(t) = ^Cx('0(O))e (70) 

X 

for some constants c^{ip(0)) depending only on the initial 
condition ■0(0). Assume now that for this initial condi¬ 
tion there is a pseudo-periodic solution with period T, 
namely, that 


■0(r) = 0(0) exp(—ia). (71) 

Comparing the eigenvector expansions of both sides of 
this equation we get the consistency condition 

(0)e-‘“ = (0)e-™>^^ for all y (72) 

which for incommensurable frequencies can be only 
satisfied if T satisfies 

T = T^:=— (73) 

wx 

for some y and simultaneously c.^{ip) = 6^^ [2l[. This 
means that for a generic matrix H, the only periodic or¬ 
bits are the ones emerging from the eigenstates of the 


2ir 


PNiE) = -^3 J dae-'<^+i) 


E: 




iTrhy/\Mpo - 1 


exp 


(74) 



Here, the pseudo-periodic orbits and their actions are still 
given by Eqns. (EH) and ED), while the stability matrix 
Mpo is given by Eq. (|T7l) but with the perpendicular co¬ 
ordinates increased by one further dimension. Also, due 
to the coincidence of the integration along the orbit and 
along the global phase, the primitive period Tppo is now 
determined by a full cycle, i.e. the smallest (non-zero) 
time for which ip{t) = 0(0). The Maslov index a can 
still be calculated according to Sec. IHI Cl 

Since the periodic orbits are isolated, our trace formula 
can be applied directly. In the non-interacting case it is 
instructive to perform the integration over a (respon¬ 
sible for selecting orbits with fixed given total number 
of particles N) exactly. The pseudo-periodic orbits are 
organized in L families corresponding to the L different 
eigenvectors of the matrix H. Consider first the primitive 
pseudo-periodic orbit associated with the eigenvector y, 
whose time-dependence can be explicitly constructed as 
in Eq. (15^ and therefore has frequency w^. As charac¬ 
teristic of harmonic systems, the period is independent 
of the energy E of the trajectory, which is given by 


E = 


E (”7’) 




Xx) 


- oE^X- 


(75) 


Note that the energy E appearing in the trace formula 
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and as the argument of the many-body DOS has noth¬ 
ing to do a priori with the eigenvalues of the single¬ 
particle problem, beyond the fact that the spectrum of H 
is part of the parameters that define the many-body prob¬ 
lem. From the point of view of semiclassics in second- 
quantized systems, in the non-interacting case the single 
particle energies simply provide the frequencies of 
the harmonic problem that defines the classical limit. 

Using Eqs. m and (1431) the action of the fcth repe¬ 
tition of any member of the y-family is easily found to 
be 




(a -I- 2'nk)h 


rix) 


(76) 


and therefore, using Eq. GU we get the action 




a 27rfc 




(77) 


and the period 


tW(e) = 


a + 2'Kk 


Wy 


(78) 


of the pseudo periodic orbits in terms of the energy. Fol¬ 
lowing the discussion below Eq. (El), the period of the 
primitive pseudo orbits is then given simply by 21 ^jw^. 

The next step is the calculation of the stability ma¬ 
trices and Maslov indexes entering the trace formula 
Eq. (l7^ . This is a standard exercise for harmonic sys¬ 
tems, and we present it here only to illustrate the concep¬ 
tual relation between eigenstates of single particle prob¬ 
lems in first quantization and the semiclassical approach 
to second-quantized many-body systems of indistinguish¬ 
able particles. 

The stability matrix is given by the local properties of 
the classical evolution around a specific pseudo periodic 
orbit, Eq. dMl), as the linear transformation relating small 
initial and final deviations from the 

reference orbit. Using the linearity of the classical dy¬ 
namics we easily get 

<5u(^)(T«) = (79) 

— ylx ) . du^^^(0)e ex K 

x' 


In view of Eq. (gZl), we then obtain for the components 
of the deviations along the directions y' = x~‘~ perpen¬ 
dicular to the y-orbit 

))e“ = (80) 


This equation indicates that the stability matrix is obvi¬ 
ously block-diagonal, with the block corresponding to y' 
being simply a rotation matrix with angle 


e 


X,(k) 

x' 


(a + 2TTk)— -a. 


and therefore 
- 1 


n 

X't^X 


det 


cos0^;^^^ — 1 
-sin0^;('=^ 


sin 0 ^;('=) 
cosd^;('=)-i 


This yiealds eventually 


Aff ^ - 1 


x't^X 


sin 


a + 27rk ey 


(82) 

(83) 


Moreover, the Maslov index for the /c-th repetition of the 
y-family is given by ( 2 ^ 


x'/x 


(k + 

-1 

^x' _ 

_V 

2ttJ 

®X 


+ 1, (84) 


where [xj denotes the integer part of x. 

Substitution of the actions, stabilities and Maslov in¬ 
dexes in Eq. El) leads to the semiclassical trace formula 
for the free bosonic field with fixed total number of par¬ 
ticles N as a sum over pseudo-periodic orbits and their 
repetitions: 



If necessary, the last integration over a can be performed 
by taking a = 0 in all smooth terms, and calculating 
exactly the integral involving the highly oscillatory ones, 
and this may be indeed the way to proceed for specific 
calculations based on the pseudo-periodic orbits. 

Eor the smooth (Weyl) contribution we left out again 
the a-integration and easily get 


P.iE)= I 


277 

J dae-'°‘^ 


0 



d^pS 


E - 


( 86 ) 


a(q^+P^) 


Eor completeness, we show below the consistency of the 
trace formula, Eq. (|85p . with the semiclassical quantiza¬ 
tion procedure for direct quantization of invariant man¬ 
ifolds in phase space, so-called EBK-quantization, valid 
only for integrable systems. In the following we apply it 
to the case of interest here. 

Within EBK quantization, the classical Hamiltonian is 
first written in terms of a new set of canonical variables 
(</>,J) where J = J(n) are combinations of the classi¬ 
cal constants of motion, such that H{4>,3) = H{3). In 
our case, these functions are simply given by J = hn(xp) 
with defined in Eq. (IMl) . In the new variables the 

Hamiltonian is given by 


^(J) = E(4-|)^x- 


(81) 


(87) 
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In a second step, EBK quantization selects the values of 
the classical actions J such that 

J^=h ( 88 ) 


In order to transform this sum over quantum numbers 
corresponding to quantized occupations into a trace for¬ 
mula where periodic (or pseudo-periodic) orbits appear, 
we proceed in a similar way as in the semiclassical quan¬ 
tization of harmonic systems and introduce 


with = 0,1, 2 ,... and with indexes given by the 
number of turning points of any classical trajectory evolv¬ 
ing in the phase space manifold defined by the set of 
quantized constants of motion, so in our case = 2 for 
all x- The EBK-quantized energies are then obtained by 


Z(/3, a) = dEe-^^ E ^ " E j -x 

= (93) 

n 


En^H 


3 = h 



giving for our case 


(89) 


Performing the sums given as geometric series, we get 


Z(/3,a) = ]^ 
x' 


(94) 


En — 

X 


(90) 


which is suitable to compute the inverse Laplace trans¬ 
form required in Eq. (1^ 


providing a proper and physical interpretation of the 
EBK quantization condition in the context of the free 
bosonic field: quantization of the many-body energy lev¬ 
els is due to quantization of the occupation numbers. 
As it is well known Q, the EBK quantization of linear 
(harmonic) systems is exact, and indeed this is the exact 
quantum mechanical spectrum of this system. 

So far, we used a well known quantization method 
to derive a well known result in the framework of first- 
quantized systems. What makes the second-quantized 
approach rather special is the status of the phase-space 
observable representing the total number of particles N. 
The key point is that the conserved quantity 

= (91) 

i 

plays a distinctive role in the field theoretic scenario, 
namely, it labels subspaces of given total number of par¬ 
ticles. Note that such a condition is never encountered 
in the description of particle (instead of field) systems. 
There, the function N{ip) is simply the sum of actions, 
and there does not exist a physical interpretation as an 
observable. 

This detail makes the semiclassical approach for fields 
with conserved number of particles conceptually differ¬ 
ent from its interpretation as a set of first-quantized har¬ 
monic oscillators. In particular, if one wants to study the 
many-body spectrum with a given, fixed one must 
project the EBK DOS. This is again accomplished by in¬ 
troducing a variable a, playing the role of an imaginary 
chemical potential, 

2tt 

Pn{E) = ^J dae-^^^ (92) 

0 

n \ X / 


( 95 ) 

n \ X / 

-| ^e—ioo 

= Tr^ dPe^^Z{P,a), 

J e—ioo 

by means of the Bromwich formula where e is real and 
positive and chosen such that it is larger than the real 
part of all the poles of Z{P,a). These poles are easily 
found to be located at 


, fc = 0,±l,±2,... (96) 

and are naturally labeled by the index x of the single¬ 
particle energy that, as we know, also denotes periodic 
orbits. The second index k labeling the positions of the 
poles is naturally interpreted within the trace formula as 
the repetition of the pseudo-periodic orbits as well. 

Having at hand our trace formula, Eq. (1851) . we are 
ready to identify and give physical interpretation to the 
different factors appearing in the “exact” trace formula 
obtained by evaluating the contributions from the residua 
of Eq. (IMl) at the poles given by Eq. (1%)) . Consider first 





Q: + 27rfc 


E 


(97) 


giving the energy-dependent term in the action of 
Eq. dHSl). Now 




0=0. 


(k) 


1 _ 


1 -1 


2 e-y 2 , 


2 isin(ii^^ - f 


(98) 


give the stability factors in Eq. (1551) . Putting all together. 
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we finally get 


(99) 


1 “ 

P 


X ^ fc 


=-oo nx'#x2isin( 


Q:+27rfc ^x' 

2 ev 


)' 


in full agreement with our trace formula, Eq. (IMI) . Note, 
however, that the terms with k = 0 must be computed 
independently, and in the limit a = 0, appropriate for 
the asymptotic regime iV —>■ oo, they precisely provide 
the Weyl term, Eq. (|55)) . 


V. EXTENSIONS AND RELATION WITH 
PREVIOUS APPROACHES 

The semiclassical trace formula provides a fundamental 
and rigorous connection between the spectrum of a quan¬ 
tum system and the features and properties of its classi¬ 
cal limit. As such, it should be able to cover less general 
approximations that rely on classical information to ex¬ 
plain features of the quantum mechanical DOS. Note, 
however, that the trace formula does not associate di¬ 
rectly dynamical properties of the classical system with 
individual energy levels of the quantum spectrum. The 
correct association is between periodic solutions of the 
classical mean field equations and Fourier components of 
the full many-body quantum DOS as a function of the 
energy. Semiclassically, the emergence of discrete ener¬ 
gies is an interference phenomenon due to the coherent 
superposition of all these harmonics. 

In this paper, we have addressed the situation where 
the classical mean field dynamics, understood as a dy¬ 
namical system, is such that all the periodic solutions 
are isolated. Remarkably, while this is usually the case 
that only for fully chaotic dynamics, the important case 
of a non-interacting bosonic field falls into this category 
as well. This is because the periodic solutions of the cor¬ 
responding classical limit, which is harmonic, are again 
isolated. 

As it is well known [mi a generic dynamical system 
is actually neither integrable nor chaotic, and the trace 
formula and more generally, semiclassical quantization 
methods for hrst-quantized systems were correspondingly 
generalized in order to describe also the integrable-to- 
chaotic transition [11. Also, extensions and generaliza¬ 
tions of the trace formula can be used to quantize se¬ 
lected, specific locally integrable dynamics of the clas¬ 
sical phase space. In the following we present a brief 
discussion of the connections between the approach pre¬ 
sented here and other semiclassical methods aiming to 
associate features of the many-body quantum spectrum 
with special classical structures. 

The first, non-generic, situation appears if the classical 
mean-field equations admit a static solution defined 

by 


ifi0(f) 






drp*{t) 




= 0 


( 100 ) 


such that the classical mean-field motion is strictly stable 
in all directions around In this case, the spectrum 

of the system contains energy levels associated with the 
quantization of the normal modes describing small oscil¬ 
lations around the stable fixed point. Long-wavelength 
oscillations of the DOS due to this sequence of energy 
levels will be described by our trace formula for the har¬ 
monic system obtained by the quadratic expansion of 
(0*(0^0(0) around 0(0 = In the con¬ 

text of Bose-Hubbard models, this result corresponds to 
the well known Bogoliubov approximation (^ . and by 
its very construction it is valid only for energies E such 


that 


~ E. 


Generic dynamical systems display a mixed phase 
space with a mixture of regular (locally integrable) and 
chaotic motion, often with complicated fractal structures 
along their borders [l^. None of the versions of the trace 
formula presented here are strictly valid for this situation. 
However, the whole machinery developed to incorporate 
bifurcation effects characteristic of the mixed dynamics 
scenario in the trace formulas for first-quantized systems 
Q can be directly imported into the many-body context. 
The semiclassical quantization of regular islands in phase 
space can be performed equivalently by constructing local 
constants of motion and applying the EBK quantization 
locally [T^ . These approaches may be very important in 
the field theoretical context, as numerical investigations 
consistently indicate that the phase space of the dynam¬ 
ics induced by the mean field Hamiltonian has a very 
complicated structure, where hard chaos might be very 
difficult to observe 


A special situation arises for the case L — 2 where in¬ 
teractions are present, rendering the classical dynamics 
nonlinear, but the system is still integrable . Here exist 
two constants of motion (the total energy and number 
of particles), hence as many as degrees of freedom. This 
two-site Bose Hubbard model has been extensively stud¬ 
ied both quantum mechanically and classically, and can 
be exactly mapped into the Josephson Hamiltonian de¬ 
scribing bosonic excitations in superconducting devices 
[^ . Semiclassically, this system has been extensively 
studied by means of the EBK quantization method ap¬ 
propriate for classically integrable systems. The con¬ 
struction of the classical actions and the quantization 
conditions is, however, substantially more complicated 
[ttI - I^ . In order to apply the methods presented here 
to this situation, it is convenient to work directly in 
the reduced phase space obtained by fixing the total 
number of particles. The dynamics is now essentially 
one-dimensional, and our trace formula can be used, as 
energy-conserving one-dimensional systems possess both 
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integrable dynamics and isolated periodic orbits. The re¬ 
sult of this calculation is just the WKB approximation 
to the energy levels, which can be improved in several 
ways [281 l30l|. This equivalence between trace formulas, 
EBK quantization and WKB methods for conservative 
one-dimensional systems is well known in the context of 
first-quantized semiclassics Q. 

Furthermore, our trace formula allows one to calcu¬ 
late the contribution to the DOS from periodic orbits 
with both stable and unstable local classical cKnamics, 
a generic case for multi-dimensional systems |^. Note 
that for the non-interacting limit of the Bose-Hubbard 
model, we encountered already the situation where the 
local classical flow around the isolated periodic orbits is 
stable. The difference between the character of the differ¬ 
ent degrees of freedom is fully encoded in the properties 
of the stability matrices around the classical periodic or¬ 
bits, and therefore the application of the trace formula 
requires their explicit calculation. The spirit of this ap¬ 
proach is not the complete enumeration of all periodic 
orbits, but the study of the contribution of particular so¬ 
lutions to the DOS. The study of such effects has been 
successfully carried out in the first-quantized approach 
to the helium atom and to the semiclassical quanti¬ 
zation of solitons in the context of continuous quantum 
field theory 321. 


VI. CONCLUSIONS 

We have presented a rigorous derivation of the semi¬ 
classical approximation for the quantum mechanical DOS 
of many-body quantum systems described by bosonic 
quantum fields on finite lattices, starting from the exact 
path integral form of the many-body propagator. We 
showed explicitly how to derive both the smooth (Weyl) 
and oscillatory (Gutzwiller) contributions to the DOS, 
and provide a trace formula for the later. Our approach 
follows and generalizes the original pioneering methods 
introduced by Gutzwiller for chaotic single-particle sys¬ 
tems. We avoided the coherent state representation, with 
its characteristic need to complexify the classical limit of 
the theory, by using quadrature states of the field. As 
a special feature of the field scenario, the classical limit 
is a discrete classical field and its isolated (mean-field) 
periodic solutions are the input of the trace formula. 

Another special property is the existence of a contin¬ 
uous symmetry related to the conservation of the total 
number of particles in closed systems. We applied the 
methods of symmetry-projected semiclassical densities 
of states to get an expression for the many-body DOS 
within each sector with fixed total number of particles. 
Interestingly, due to the fact that the quantum problem 
in the non-interacting case reduces to a harmonic system, 
our trace formula is applicable since the periodic orbits 
can be shown to be isolated in the generic case where the 
single-particle energies are not commensurable, as in the 
chaotic case. 


As for the single-particle case, our trace formula shows 
how the existence of discrete many-body energy levels 
emerging from the continuous, smooth background given 
by the Weyl law, is an interference phenomenon. The cor¬ 
rect density of states and its characteristic profile made 
up from Dirac-delta peaks is built by the coherent effect 
of oscillatory contributions, one for each periodic orbit, 
i.e, periodic mean-field solutions. 

The study of the quantum manifestations of classical 
solitons (particular solutions of the non-linear equations) 
is matter of recent interest in the cold-atoms community. 
We expect that the application of our methods may help 
to understand and quantify interference between different 
solitonic contributions to the many-body DOS in terms 
of harmonics building up the trace formula, a pure quan¬ 
tum effect not to be confused with the possibility of hav¬ 
ing remnants of wave interference at the mean-field level. 
Specifically, if the classical non-linear field equations ad¬ 
mit two solitonic solutions with actions Si and S 2 which 
are stable or only mildly unstable, we can approximate 
(up to higher repetitions) Eq. (ICT) as 


52(E)^ (101) 

and the DOS will display a characteristic beating pat¬ 
terns due to the interference between the two oscillatory 
functions, reminiscent of (super-) shell effects in the nu¬ 
clear physics context 00- We see that this kind of 
analysis will provide a deeper, semiclassical understand¬ 
ing of long-range structures in the energy spectrum of 
many-body systems. 

Our work paves the way to a systematic study of the 
role of classical field solutions in the many-body DOS 
for discrete quantum fields. The close formal analogy 
of the many-body trace formula (l?T]l with the Gutzwiller 
trace formula implies, once the conservation of total num¬ 
ber of particles is taken care of, that the methods rou¬ 
tinely used in the context of single-particle semiclassics 
can be straightforwardly generalized to show that the 
spectral fluctuations have universal correlations on the 
scale of the mean level spacing if the classical field equa¬ 
tions are chaotic. The universal spectral correlations 
are linked to interference between periodic orbits with 
quasi-degenerate actions and periods beyond the Ehren- 
fest time, given by ^ log?i“^ for single-particle sys¬ 
tems with a classical limit with Lyapunov exponent A. In 
the many-body case, correspondingly, quantum interfer¬ 
ence is additionally governed by another log-time scale, 
the Eherenfest time A“^ log N. Here, A denotes the aver¬ 
age Lyapunov exponent of the assumed non-linear mean- 
field dynamics of the classical field limit. Many-body 
interference evolves at time scales beyond this Ehrenfest 
time and therefore is not accounted for in usual so-called 
Truncated Wigner approaches [s^. 

Interfering quasi-degenerate periodic mean-field solu¬ 
tions are expected to lead to the emergence of universal 
many-body spectral fluctuations. While the close formal 
similarity to the chaotic single-particle case directly im- 


Pn{E) ~ Ai cos ( -Si{E) -h A 2 cos ( - 
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plies corresponding RMT-type expressions for the spec¬ 
tral many-body correlator or form factor, interesting new 
parametric correlation functions, involving (changes in) 
particle number an interaction strength, appear for the 
many-body case. 

The non-interacting case has, however, non-generic 
spectral fluctuations that do not correspond to the ex¬ 
pected Poissonian spectra of integrable systems, a pecu¬ 
liar consequence of the field theoretical scenario where 
the free field is actually a linear, not only integrable, 
system. Remarkably, although for non-interacting sys¬ 
tems the linearity of the classical limit makes it highly 
non-generic from the point of view of dynamical sys¬ 
tems, this is the generic case for free fields from the 
field-theoretical side. Exploring the application of semi- 
classical techniques for the quasi-integrable case where a 
small interaction is treated within semiclassical pertur¬ 
bation theory demands then techniques adequate to this 
non-generic situation [s^. 

With our work we hope to contribute to the qualitative 
and quantitative analysis of both universal and system- 
specific features in the energy spectra of many-body sys¬ 
tems by means of periodic orbit theory. 
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Note: A derivation of a related trace formula has 
been performed independently by R. Dubertrand and 
S. Muller using a periodic orbit theory approach in 
the reduced phase space, see Sec. El 


Appendix A: Derivation of the semiclassical 
prefactor of the trace formula 

The simplification of the product of the two determi¬ 
nants appearing in Eq. (1441) to the form in which they 
appear in Eq. c an, to a large extent, be performed 
by following Ref, [l^ (see Eqs. (2.14-3.12) there). There¬ 
fore, here only the main steps will be presented, in order 
to show how to correctly account for the phase difference 
a. Eor notational simplicity, in the following we will use 
t instead of to represent the period of the pseudo- 
periodic orbit 7 . 

Eirst, guided by the Gutzwiller trace formula, where 
the semiclassical prefactor is determined by the mon- 
odromy matrix, we try to bring the second matrix in 
Eq. (1441) into a form, where ■j/j appears in such a form 
that without derivatives it could be replaced by p and q. 
This is achieved by 


one in order to obtain (compare with Eq. (3.4) therein) 


/ I sin a 9Qi/’j^(0) _ 1 

I 5q_L ~ COSO; Op COSO; 

I dR'tp{t) _ 1 dR'tp{t) I sin a. 

\ Oq L cos a dp ~ cos a. 




(cos a) 


X det 


-(2L-2) 


det 


'd (5R^||(t)e-*“,E, 


d ( 0 ) 

d - p^, jR^^(t)e~‘° - qu, 6 >) 

d{q±,^il’± (0), A^) 


(A2) 

where 9 is the initial global phase of the nonlinear wave 
and N-y is the number of particles determined by the 
nonlinear wave. 

Next, in the first of these two matrices, we again re¬ 
place the derivative with respect to ( 0 ) by one with 
respect to S {ip±{0) exp (ia)), 


det 


(9(qj.,Si/)_L (0), A..^) 


, .L-2 , , '9(3^u(0)e*“-pj_,Ki/>^(t)e *“-qu,6') 

(cos a) det- 


( 0 )ei“, A..y) 


(A3) 


Thus, 


I sin a OSfi/’j.CO) 1 

rlof I COSO Op COSO 

I ddllPit) _^ I sin a 

\ cos ot Op ~ cos o 




= (cos a) 


X det 


-L 



1 det 1 

{ 95R| 

(^11 We 

1 9(A,A^)\ 



d{t,0) 

93i/)||(0) 1 


d (jR'0^(t)e - q_L, pu - Si/>_L(0)e‘“) 


9(q±,3’0_L ( 0 ) e'“) 


= — (cosa)~ 


46 ^ / 36 


dN^ 




det 


(^11 We 

d{t,0) 


X det 


3 (K'0_L(t)e - q_L, pj_ - 3'i/)_L(0)e‘“) 

3{q±,Qipj_ (0)ei“) 


(A4) 


The additional minus sign in the last step results from 
the fact that 


det 


. sing d^'4’±{0) _ 1 \ 

Oq_L ~ cos Ol Op cos o I _ 

OjRi/>(t) _ 1 I sin OL j 

Oq_L cos a Op cos a / 


det 


1 i9 (5'j/)_L(0)e‘“ — pu, Kt/’(t)e — q) 
cos a 5(q±)P) 


(Al) 


With this the steps in Ref. [l9| can be carried over one to 


3^ (i/>ll (t)e-‘“) 
90 


5^11 ( 0 ) 


1 / 9A^ \ 

2 [dQ^„(0)l 


(A5) 

The remaining prefactors stem from the special form of 
the equations of motion, Eq. (EH), together with the re¬ 
lation between rp and q, p, Eq. (H51) . 































Appendix B: Uniqueness of the phase difference of 
psendo-periodic orbits 

In this appendix, we show that if 7 is quasi-periodic 
with phase difference a, i.e. 

rpiT^) = ■0(0) exp(-ia), (Bl) 

with pseudo-period 

Tppo = T^/m, m G N, (B2) 

then there is no /3, which is not an integer multiple of 
a/m, such that for some time t* the trajectory satisfies 

0(r*)=0(O)exp(-i;9). (B3) 

To see this, assume there would be such a phase /3 and 
time T* < Tj. Since P is not an integer multiple of a/m, 
T* also is not an integer multiple of Tppo Then there is 
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also an integer 

m*T* < 0 < (to* -f l)r* (B4) 

and the trajectory would have to satisfy 

0(Ty) = 0 (to*T*) exp(—i(a — m*P)), (B5) 

i.e. the trajectory would be pseudo-periodic with phase 
difference a — m*P and primitive period Tj — m*T* < T*, 
which is not an integer multiple of Tppo- Now we can 
replace /3 by a — m*P and repeat the argumentation again 
yielding an even smaller primitive period Ty — m*T*. 
Since to is supposed to be the largest possible number 
such that 

xf){Tj/m) = 0(0) exp(—ia/TO.) (B6) 

this procedure can be repeated infinitely often yielding 
finally T* = 0. However, having a pseudo-period T* = 0 
means, that 0(t) is an eigenstate of the mean-field Hamil¬ 
tonian, and therefore an orbit of zero length, which are 
not included in the oscillatory part. In the same way, 
one can show that the case T* > Ty leads to the same 
contradiction by simply replacing the roles of a and /3. 
Thus, Eq. (11511 is indeed the correct result. 
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